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Erefacs 

The  purpose  of  this  study  was  to  develop  a  computer  algorithm  capable  of 
generating  correlated  random  processes  for  simulating  time-varying  samples  of  chaff  cloud 
radar  cross  sections.  1'he  need  for  this  model  arises  in  realistic,  dynamic,  electronic 
warfare  training  scenarios.  The  model  may  also  be  used  to  evaluate  radar  systems  and  their 
effectiveness  against  chaff. 

The  algorithm  was  developed  on  a  VAX  11/785  computer  using  a  4.2  BSD  version 
of  the  UNIX  operating  system,  the  FORTRAN  77  programming  language  and  The 
International  Mathematical  and  Statistical  Library  (IMSL)  of  subroutines. 

I  developed  an  interest  in  electronic  warfare  and  chaff  modeling  through  some  of 
my  classes  at  AFIT  and  was  fortunate  to  find  a  topic  that  matches  my  interests  and  an 
advisor  who  was  enthusiastic  about  supportir.0  it. 

Richard  P.  Fray 
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In  this  thesis,  time-varying  radar  cross  sections  of  chaff  clouds  are  generated  for 
use  in  radar/ECM  computer  simulations,  under  the  assumption  that  scattering  from  chaff 
clouds  is  a  wide  sense  stationary  tandom  process.  For  a  jointly  gaussian  random  process, 
second  order  statistics  are  developed  from  correlated  samples  of  radar  cross  sections. 
Applicable  statistical  tests  are  performed  to  validate  a  set  of  generated  samples.  The 
goodness  of  fit  of  the  samples  compared  to  the  prespecified  density  function  is  determined 
along  with  a  comparison  of  the  correlation  coefficients  from  the  generated  set  to  the  desired 
correlation  coefficients.  Rayleigh  and  Weibull  statistics  are  also  derived  from  gaussian 
variables  to  present  an  alternative  way  of  describing  the  probability  distributions  of  chaff 
cloud  cross  sections.  Topics  for  further  study  are  suggested. 


SIMULATION  OF  CHAFF 
CLOUD  SIGNATURE 


Simulations  of  radar,  electronic  countermeasure  (ECM),  or  electronic 
countercountermeasure  (ECCM)  systems  involving  chaff  have  been  of  interest  to  the 
military  community  for  years  in  the  evaluation  of  such  systems  and  in  electronic  warfare 
scenarios  used  in  training  pilots,  radar  operators  and  electronic  warfare  officers. 

The  term  chaff  denotes  a  confusion  type  of  radar  ECM  which  employs  a  large 
number  (hundreds  of  thousands)  of  resonant  dipoles  packaged  into  small  cartridges  (16:1). 
The  individual  chaff  dipoles  are  manufactured  in  a  variety  of  ways  and  from  different 
materials  including  glass  fibers  dipped  in  molten  silver  or  aluminum,  but  the  standard 
material  used  in  the  production  of  chaff  is  aluminum  foil  cut  into  thin  strips  (5:18-14). 

Individual  chaff  dipoles  have  bandwidths  corresponding  to  about  10%  of  their 
resonant  frequency,  so  a  single  cartridge  of  chaff  usually  contains  several  different  lengths 
that  will  resonate  over  a  wide  band  of  frequencies.  When  the  tubes  are  dispensed  into  the 
air  from  an  aircraft,  helicopter,  or  as  a  mortar  round,  the  cartridges  burst  and  create  a  "radar 
smoke  screen"  that  confuses  the  enemy  radar  (5:18- 16). 

The  radar  echo  produced  is  the  result  of  the  return  from  a  large  number  of 
individual  dipoles  (5:18-10).  The  echo  can  be  used  to  saturate  a  radar  coverage  area,  break 
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lock  on  a  tracking  radar,  or  create  false  targets  (5:18-17).  Corridors  of  chaff  can  also  be 
formed  to  provide  friendly  aircraft  with  areas  to  fly  through  without  being  detected  or 
tracked  by  enemy  radar. 

Background 

Various  attempts  have  been  made  to  simulate  chaff  clouds  realistically  using  the 
various  parameters  associated  with  the  clouds  and  their  probability  distributions.  Studies 
have  been  made  to  analyze  the  distributions  of  dipole  angle  orientation  (21),  radar  cross 
section  (1,  2),  and  drift  velocity  (3).  The  existing  methods  of  synthetically  generating  radar 
reflections  from  chaff  clouds  are  either  deterministic  or  at  best  restricted  to  the  first  order 
statistics  which  merely  indicate  the  magnitudes  of  the  fluctuations  of  the  cross  sections 
about  some  mean  value.  The  deterministic  methods  predict  the  exact  outcome  in  any  given 
situation  and  often  require  some  very  involved  calculations.  Due  to  limited  computer 
memory  and  time,  only  a  few  hundred  dipoles  can  usually  be  analyzed  in  representing  one 
chaff  cloud  (14:1). 

The  restriction  to  first  order  statistics  is  not  very  satisfactory  since  in  a  dynamic 
situation  one  has  to  know  not  only  the  magnitudes  but  also  the  rates  of  the  fluctuations.  On 
die  other  hand,  statistical  methods  are  >t  concerned  with  any  particular  situation.  They 
predict  by  a  probabilistic  method  what  might  happen  under  a  given  set  of  conditions.  There 
is  no  restriction  on  the  number  of  dipoles  a  cloud  can  consist  of.  In  fact,  the  larger  the 
number,  the  more  accurate  the  representation  will  be  (14:1). 

To  determine  the  rates  of  the  fluctuations,  second  order  statistics  of  the  chaff  cloud 
radar  cross  sections  are  needed.  These  statistics  give  the  probability  of  jointly  finding  two 
values  of  a  random  variable  at  different  times. 
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Tlic  radar  return  of  a  chaff  cloud  at  one  instant  of  time  will  be  correlated  to  some 
extent  with  the  return  at  a  previous  time.  For  small  time  differences,  the  correlation  will  be 
high,  and  for  large  ones,  the  correlation  will  be  low  (5:18-8).  To  model  this,  the 
correlation  function  will  be  specified  as  an  exponential  function  of  the  time  difference.  The 
purpose  of  this  thesis  is  to  develop  an  algorithm  to  generate  pseudorandom  samples  of 
radar  cross  sections,  with  these  correlation  coefficients,  for  a  representative  chaff  cloud. 

For  a  given  autocorrelation  function,  second  order  stadstics  will  be  generated. 

This  model  will  remain  a  fairly  simplified  one  in  the  respect  that  it  will  deal  only 
with  the  radar  cross  sections  of  a  chaff  cloud.  The  aerodynamic  and  electromagnetic 
qualities,  drift  velocities,  and  birdnesting  effects  (clumping  together  of  several  chaff 
elements)  of  individual  dipoles,  and  the  transient  effects  of  cloud  blooming  immediately 
following  dispensation  of  a  chaff  cartridge  will  not  be  addressed. 

Assumptions 

The  primary  assumption  made  is  that  radar  backscattering  from  chaff  clouds  is  a 
wide  sense  stationary  random  process.  This  implies  that  the  statistics  of  one  dipole  or 
cloud  will  not  change  from  one  time  instant  to  another,  and  that  the  autocorrelation  is 
dependent  only  on  the  time  difference  at  which  the  dipole  or  cloud  is  evaluated.  It  also 
implies  that  the  average  or  mean  value  of  the  statistics  is  independent  of  time.  This  is 
justified  due  to  the  fact  that  chaff  has  a  very  high  aerodynamic  drag  and  therefore  will  not 
alter  its  position  very  much  between  radar  returns. 

Each  elemental  dipole  will  be  assumed  to  have  random  phase  and  amplitude  for  two 
reasons.  The  orientation  of  the  dipole  may  change  due  to  rotation  of  the  dipole,  and  the 
distance  between  the  radar  and  the  dipole  center  may  change.  This  implies  the  cloud  will 
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also  have  random  phase  and  amplitude.  Finally,  it  is  assumed  the  separation  between 
dipoles  is  great  enough  to  prevent  any  mutual  coupling  or  interaction  between  adjacent 
dipoles.  The  effect  of  the  random  phases  and  amplitudes  merely  suggests  that  the  model 
will  generate  a  radar  cross  section  that  is  uniform  throughout  a  chaff  cloud. 

General  Approach 

The  model  entails  generating  a  set  of  correlated  psendorandom  samples  from  a 
gaussian  probability  density  function.  Two  exponential  correlation  functions,  as  shown  in 
equations  1  and  2,  will  be  specified  for  each  set  of  samples. 

P(t)  =  exp  [  -  |t|  3  (1) 

p  (t)  =  exp  (  -  t2  )  (2) 

Where  t  corresponds  to  the  pulse  repetition  interval  of  a  radar  system  and 
represents  the  time  between  consecutive  samples  of  the  radar  cross  section.  Several  sets  of 
data  will  be  generated  with  different  quantities  of  cloud  samples,  for  both  correlation 
functions,  to  form  a  basis  of  comparison  and  validation  for  the  model.  The  number  of 
samples  will  be  limited  to  at  most  100  because  in  using  pulse  integration  to  obtain 
additional  data,  rarely  arc  more  than  100  pulses  integrated. 

For  each  set  of  samples,  a  chi-squated  goodness  of  fit  test  will  be  performed  to  see 
if  the  generated  samples  do  indeed  fit  the  theoretical  distribution  they  were  originally 
generated  from.  Then  a  comparison  will  be  made  between  the  autocorrelation  coefficients 
of  the  generated  samples  and  the  originally  specified  exponential  correlation  functions. 
Finally,  Rayleigh  and  Weibull  samples  will  be  derived  from  the  gaussian  samples. 


First  the  existing  theories  of  representing  chaff  cloud  returns,  methods  developed 
for  generating  correlated  random  variables,  and  the  derivations  of  Rayleigh  and  Weibull 
statistics  from  gaussian  variables  are  discussed  in  chapter  2. 
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Then  the  algorithms  for  generating  the  correlated  random  samples  is  flowcharted, 
developed,  and  discussed  in  chapter  3. 

Chapter  4  presents  validation  of  the  model  through  applicable  statistical  testing  of 
the  generated  sets  of  data. 

Finally,  conclusions  are  made  and  recommendations  for  further  study  suggested  in 
chapter  5. 

A  review  of  the  applicable  probability  theory  and  a  computer  listing  of  the  algorithm 
are  presented  as  appendices. 
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II.  Existing  Theory 


In  this  chapter,  existing  theories  developed  for  generating  correlated  random 
variables  and  some  metiiods  used  to  represent  chaff  clouds  are  reviewed.  Also  the 
derivations  of  Rayleigh  and  Weibull  statistics  from  gaussian  variables  are  shown. 


Chaff  Models 

One  model  reported  in  the  literature  was  developed  at  the  Georgia  Institute  of 
Technology  in  1977  (4).  A  preliminary  chaff  radar  cross  section  (RCS)  model,  developed 
to  estimate  the  RCS  of  deployed  chaff,  was  unproved  and  extended  to  compute  the  mean 
RCS  of  deployed  chaff  bundles  and  their  variations  about  the  mean  value  known  as 
scintillations  (4:1).  The  Georgia  Tech  study  proposed  two  methods  to  simulate  the 


problem  of  chaff  pattern  backscaiter  scintillation,  the  Rayleigh  distributed  draw  meuuxi, 
and  the  dynamic  time  series  model  (4:3).  Each  model’s  major  features  will  be  discussed 
briefly. 


Backscatter  power  from  a  chaf  f  cloud  is  not  deterministic.  Rather,  it  varies  around 
some  mean  value.  It  can  be  thought  of  as  an  analog  process  being  digitally  sampled  by  a 
pulsed  radar  system  (4:10).  The  Rayleigh  draw  model  entailed  inputting  parameters  such 
as  the  RCS  mean,  the  period  of  the  highest  scintillation  rate,  and  the  observation  time 
desired.  A  draw  was  then  made  from  a  Rayleigh  distributed  random  number  generator  and 
processed  through  an  interpolation  formula.  The  interpolation  formula  allowed  the  user  to 
determine  the  RCS  at  any  instant  of  time  between  samples  or  at  any  rate  slower  than  the 
sampling  rate.  The  resulting  output  was  an  RCS  chaff  pattern  generated  as  a  function  of 


time. 

Amplitude  fluctuations  have  been  shown  to  fit  a  power  Rayleigh  distribution  (15), 
so  this  method  is  realistic  for  most  chaff  simulations.  Other  distributions  may  be  used  in 


different  applications  with  satisfactory  results.  To  ensure  accurate  sampling,  draws  from  j 
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the  Rayleigh  generator  must  be  made  at  the  Nyquist  rate  or  higher,  corresponding  to  the  j 

highest  scintillation  rate.  But  care  must  be  taken  because  inputting  a  rate  higher  than  the  I 

maximum  scintillation  rate  would  introduce  errors  due  to  the  inconsistency  of  specifying  a 
maximum  and  then  exceeding  it 

A  problem  may  also  arise  with  the  validity  of  the  generator  draws  if  the  statistics  of 
the  RCS  change  with  time.  During  the  bloom  period  this  may  present  a  serious  problem, 
but  for  first  order  models,  if  the  mean  RCS  was  assumed  to  vary  slowly  with  respect  to  the 
rate  at  which  a  new  RCS  value  is  generated,  then  the  results  would  be  considered  valid 
(4:14-18).  No  testing  was  done  to  compute  the  correlation  of  the  RCS  from  pulse  to  pulse, 
so  no  comparison  could  be  made  between  the  correlation  of  the  generated  variables  and  the 
desired  correlation. 

A  more  mathematically  based  method  developed  in  this  Georgia  Tech  report  was  the 
dynamic  time  series  model.  A  time  series  is  a  sequence  of  successive  observations  or 
simulations  of  a  stochastic  or  random  process.  The  observations  are  usually  equally  spaced 
and  could  represent  the  RCS  return  of  a  chaff  cloud  from  consecutive  radar  pulse  returns 
(4:19).  The  time  series  models  basically  use  previous  samples  to  determine  future  samples. 

This  model  is  more  rigorously  developed  but  could  require  longer  computer  processing 
time  than  the  Rayleigh  draw  method  (4:34). 

Another  chaff  model,  also  developed  at  Georgia  Tech,  was  the  Chaff  Theoretical/ 

Analytical  Characterization  and  Validation  Program  (7).  This  model  incorporated  the 
aerodynamic  and  electromagnetic  behavior  of  chaff  dipoles  and  modeled  the  propagation  of 
received  signals  through  a  radar  system  (7:1). 

This  model  was  a  very  detailed  one,  taking  into  consideration  dipoles  of  different 
lengths,  and  weighting  chaff  returns  according  to  a  dipole's  angular  location  within  the 
antenna  radiation  pattern  (7:1-2).  The  approach  was  to  develop  a  quasi-deteiministic 
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description  of  the  cloud.  Further  asstunptions  made  concerning  the  dipole  density  function 
resulted  in  integrals  that  could  not  be  evaluated,  not  even  numerically.  These  integrations 
were  too  time  consuming  on  a  computer  and  that  made  the  model  unworkable. 

They  suggested  future  efforts  be  devoted  to  making  approximations  of  the  volume 
integral  and  empirical  descriptions  of  the  dipole  density  variations  (7:4).  Also  better 
experimental  data  was  needed  to  help  predict  and  validate  the  model. 

Other  studies  exist  taking  different  approaches  including  Wickliff  s  and  Garbacz’ 
(20),  which  included  mutual  coupling  among  dipoles,  and  Schiffs  (18),  which  modeled 
chaff  decoy  placement  in  relation  to  antiship  missile  defense  problems  of  the  Navy. 
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Several  methods  were  reported  in  the  literature  for  generating  correlated  random 
number  sequences.  Some  of  the  algorithms  were  strictly  designed  to  generate  gaussian 
variables  while  others  could  be  used  to  generate  any  probability  distribution.  One  method, 
developed  by  Nawathe  and  Rao,  used  the  theory  of  optimal  linear  prediction  to  generate 
cross  correlated  random  -lumber  sequences  (11). 

The  technique  first  required  generation  of  two  random  number  sequences,  for 
example  X  and  Z,  each  having  the  desired  probability  distribution.  The  vector  Z  was  used 
as  a  parent  population  for  generating  a  correlated  vector  X  Then  for  each  element  of  the 

VOptni*  Y  v(i\  f  Y  —  ( v  (  1  \  vfn  tl  \  n  nnl  Aulntt  tl»i>  U  .in  f  1 « r*  /*«>••  4 « *■*♦■*-%•*  f  X/-  mor. 

V4A  —  L^vA/  •••  Jy  a  loivuiauuu  vi  uiv  uv/oi  nuwui  vi  x  waa 

made.  For  every  x(i),  a  search  was  made  through  Z  to  find  a  sample  which  was  close  to 
x(i),  closeness  defined  as  a  specified  mean  square  error  interval.  If  no  match  was  found,  a 
sample  was  picked  from  Z  which  had  the  lowest  mean  square  error,  away  from  x(i),  and 
the  algorithm  assigned  this  sample  to  the  corresponding  y(i)  (1 1:98). 

This  method  is  simple  and  straight  forward  and  requires  only  knowledge  of  the 
marginal  densities  of  X  and  Z •  Using  this  method,  several  sets  of  data  were  tested  with 
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different  desired  correlation  coefficients  and  probability  densities.  The  estimate  of  the 
coefficient  of  correlation  between  X  and  X  was  defined  (1 1 :97)  as 

nr 

p  =  (1/nr)  ^  {  [  x(j)  -  ]  [  y(j)  -  Py  ]  /ox  oy  }  (3) 

j=I 

A  table  of  the  calculated  correlation  coefficients  showed  the  difference  from  the  desired 
coefficients  increased  as  p  approached  zero,  and  was  on  the  order  of  only  1  to  2  percent 
for  a  correlation  of  p=.8  (1 1:99).  This  confirmed  the  statement  in  the  paper  that  results 
tended  to  improve  with  increasing  magnitude  of  the  correlation  coefficient  and  sample  size. 
The  model  was  further  validated  through  derivations  of  the  expressions  for  the  expectation 
and  variance  of  the  correlation  coefficients  (i  1:96). 

This  technique  proved  satisfactory  for  three  typical  distributions,  normal, 
exponential,  and  uniform,  and  thus  appears  applicable  to  any  continuous  distribution. 
Although  it  is  only  an  approximation  for  nonnormal  distributions,  the  results  of  the  Monte 
Carlo  simulation  here  are  valid  (1 1:101). 

Another  method,  developed  by  Li  and  Hammond,  involved  transforming  a  vector 
of  independent  gaussian  random  variables  into  a  vector  of  gaussian  variables  with  a 
specified  predistorted  correlation  coefficient  matrix  (8).  First,  a  vector  of  independent 
gaussian  random  variables  was  generated.  Then  this  vector  was  multiplied  by  a  linear 
irinsiormation  matrix  which  resulted  in  correlated  gaussian  variables  possessing  the 
predistorted  correlation  coefficients.  Finally,  a  nonlinear  transformation  was  performed  to 
create  a  vector  of  correlated  random  variables  with  the  desired  univariate  probability  density 
and  correlation  coefficients. 

The  difficulty  in  this  technique  lies  in  finding  the  nonlinear  transformation.  For  the 
desired  >robability  density  function,  a  set  of  transformations  involving  inverse  cumulative 
density  functions  and  the  error  function  are  required  (8:558).  The  calculation  of  the 
predistorted  correlation  coefficients  also  requires  complicated  integral  evaluations. 
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However,  data  produced  by  the  algorithm  passed  all  statistical  tests  and  generated 
correlation  coefficients  with  less  than  .1%  error  as  compared  to  the  desired  coefficients 
(8:559).  Therefore,  this  method  proved  itself  mathematically  sound  though  it  may  be  too 
complex  for  the  task  at  hand. 

The  technique  chosen  for  the  chaff  cloud  signature  model  was  developed  by 
Scheuer  and  Stoller  (17)  and  a  similar  study  done  by  Moonan  (10).  This  method  proved 
the  most  straight  forward  one  found  in  the  literature  and  was  well  suited  for  machine 
computation. 

The  technique  first  required  the  generation  of  a  vector  of  independent  normal 
random  variables  with  zero  mean  and  variance  of  one.  The  resulting  vector  will  also  have  a 
zero  mean  and  variance-covariance  matrix  specified  by  the  user  (17:278).  To  obtain  the 
correlated  vector,  the  independent  vector  is  multiplied  by  a  linear  transformation  matrix, 
created  from  the  covariance  matrix.  The  formulas  used  to  generate  the  transformation 
matrix  and  the  correlated  variables,  will  be  presented  in  chapter  3  along  with  a  detailed 
explanation  of  the  entire  algorithm  (17:279). 

Ravleigh  and  ibull  Derivations 

The  purpose  of  this  section  is  merely  to  present  a  method  of  obtaining  Rayleigh  and 
Weibull  distributed  random  variables  from  independent  gaussian  variables.  These  densities 
are  often  used  to  define  the  distributions  of  chaff  cloud  radar  returns,  and  are  therefore 
pertinent  to  this  study. 

First,  two  independent  identically  distributed  gaussian  random  variables,  X  and  Y, 
each  with  zero  mean,  and  identical  variances  s,  are  generated.  Then,  X  and  X.  are  each 
squared  and  added  together  to  form  the  variable  Z- 


Z  =  X2  +  X2 


(4) 


Finally,  the  variable  Z  is  raised  to  the  power  k. 
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Now,  by  specifying  k,  the  variable  W  can  be  transformed  into  one  of  several  different 
distributions.  For  example,  if  k=l,  W  becomes  exponential,  if  k=l/2,  W  becomes 
Rayleigh,  and  if  k=l/3,  W  becomes  Weibull  although  the  exponential  and  Rayleigh 
distributions  are  both  just  special  cases  of  the  Weibull  distribution.  The  derivation  of  these 
results  follows. 

The  moment  generating  function  was  used  to  transform  the  gaussian  variables  and 
is  shown  in  equation  6  for  Z=X2  (13:1 15). 


M*(t)  =  Efexp(tz)] 

=  E[exp(tx2 )] 


f(x)  exp(tx2 


)  dx 


-oo 


(6) 


After  making  some  substitutions  and  evaluating  the  expected  value  integral,  the 
result  is  simplified  to  equation  7. 


Mz(t)  =  [1/(1  -  2st )  ]1/2 


(7) 


This  moment  generating  function  represents  a  gamma  function  with  parameters 
a=l/2s  and  r=l/2  (9:195).  If  the  substitution  of  s=l  is  made,  then  Z  has  a  chi-squared 
distribution  with  one  degree  of  freedom.  This  simplification  can  be  justified  by  assuming 
the  generated  gaussian  variables  have  a  variance  of  one,  which  was  intended  when 
generating  the  gaussian  variables.  Since  X  and  Y  are  identically  distributed,  can 
also  be  represented  by  a  moment  generating  function  of  a  chi-squared  distribution  with  one 
degree  of  freedom.  Meyer  proved  that  a  sum  of  n  chi-squared  distributed  random  variables 
each  with  one  degree  of  freedom  will  result  in  another  chi-squared  variable  with  n  degrees 


of  freedom  (9:201).  It  can  further  be  shown  that  this  distribution  is  exponential  with 
parameter  a=l/2. 

fz  (Z)  =  (1/2)  exp(-z/2)  (8) 

Finally,  by  transforming  Z=X2+X2  into  W=Zk  using  the  Jacobian  (13:95),  equation  9 
results. 


fw(w)  =  f7(z)  |dz/dw| 

=  4(wl/k )  |dz/dw| 

=  (1/2)  exp[-(w1/k  )/2]  |  (1/k)  wht-kyk]  | 

=  ( 1  /2k)  wlO-k)  '*]  exp[-(w1/k  )/2],  w,  k  >  0  (9) 

The  Weibull  density  has  the  general  expression 

fw(w)  =  abwb4  exp(-awb  )  (10) 

If  the  parameters  are  set  to  b—  1/k  and  a=l/2,  equation  9  results.  Then  by  substituting 
different  values  in  for  k,  tire  densities  mentioned  previously  will  result. 
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In  this  chapter,  the  algorithm  wili  be  discussed  in  detail  to  help  the  reader  adapt  it  to 
his  own  application.  The  algorithm  was  written  in  the  Fortran  77  programming  language 
on  a  V  AX  11/785  computer.  A  flowchart  of  the  algorithm  is  shown  in  figure  1 . 

First,  several  parameters  have  to  be  specified  according  to  the  user's  desires.  This 
model  assumed  a  generic  pulsed  radar  system  for  receiving  the  cross  section  returns,  so  a 
pulse  repetition  frequency  (PRF)  was  specified,  t  was  used  in  the  program  to  designate 
the  time  interval  between  radar  pulses  and  is  equal  to  the  reciprocal  of  the  PRF  (19:2-3).  In 
this  model,  the  value  for  t  was  chosen  as  .001  corresponding  to  a  PRF  of  1000  Hz.  This 
is  a  reasonable  value  for  a  low  to  medium  PRF  radar  and  corresponds  to  an  unambiguous 
range  of  100  nautical  miles  (19:2).  Also,  the  number  of  consecutive  pulse  returns  desired, 
nr,  was  specified.  Three  trials  of  nr  equal  to  20, 50,  and  100  were  used  to  form  a  basis  of 
comparison. 

The  first  step  in  the  algorithm  required  generation  of  a  vector  of  independent 
gauss ian  random  numbers,  called  X  with  zero  mean  and  variance  one  (17:278). 
Throughout  this  algorithm,  several  subroutines  from  the  International  Mathematical  and 
Statistical  Library  (LMSL)  of  subroutines  were  utilized  (6).  These  subroutines  are  also 
written  in  the  Fortran  programming  language.  The  subroutine  GGNML  was  called  twice  to 
generate  two  vectors  of  length  nr  random  numbers.  One  was  used  to  generate  the 
correlated  sequence  and  the  second  was  used  in  conjunction  with  the  first  to  create  the 
Rayleigh  and  Weibull  variables  discussed  in  chapter  two. 

The  next  step  involved  selecting  the  desired  variance-covariance  matrix.  Two 
expressions  given  by  equations  1  and  2  were  used  to  represent  the  correlation  between  time 
samples.  These  exponential  functions  slowly  decrease  with  small  time  difference  and 


Figure  1.  Flowchart  of  Algorilhm 
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therefore  model  a  high  correlation  very  well. 

The  section  of  the  program  that  computes  the  correlated  vector  X  is  repeated  twice 
with  a  different  variance-covariance  matrix  each  time.  The  two  matrices  are  generated  by 
using  the  two  exponential  correlation  functions  described  in  chapter  one.  This  will  further 
provide  a  comparison  of  different  correlations  often  used  and  hopefully  help  the  user  in 
selecting  the  appropriate  form. 

A  transformation  matrix  was  needed  to  change  the  independent  vector  into  a 
correlated  sequence.  The  method  developed  by  Scheuer  and  Stoller  (17),  as  presented 
briefly  in  chapter  two,  was  used.  A  lower  diagonal  transformation  matrix  C  is  calculated 
from  the  variance  -  covariance  matrix  ££  using  the  following  formulas. 


C(i,l)  =  CC(i,l)/[CC(l,l)]1/2  1<  i  <;  nr 

(11) 

i-1 

C(i.i)  =  [CC(i.i)  -  X  t  (i.k)  11/2  1<  i  S  nr 

n?i 

\ - -  f 

k=l 

C(i,j)  =  {CC(ij)  -  Zw  [C(i,k)  *  C(j4c)  ]  }/C(j,j)  l<j<i£nr  (13) 
k=l 

C(i,j)=0  icj^nr  (14) 

Finally,  to  generate  the  correlated  vector  X,  the  transformation  matrix  £  is  multiplied  by  the 
independent  vector  X  as  shown  in  equation  15  (17:279). 


X(.i)  =  ^  C(i,j)  *  Y(j)  i  =  1,2,..., nr  (15) 

j=l 

The  next  section  of  the  program  calculated  the  autocorrelation  of  the  vector  X-  This  was 
done  to  determine  the  degree  of  correlation  between  the  samples.  The  autocorrelation  was 
computed  for  nr  time  lags  (the  entire  sequence)  througn  use  of  equation  16  (22:69). 
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i  =  1,2,. ...nr 


(16) 


nr-i 

Rx(i)  =  (1/nr)  Xx(j)X(j+D 

j=l 

The  power  spectral  density  (PSD),  one  of  the  second  order  statistics,  was  also  calculated 
from  the  correlated  sequence.  The  purpose  of  this  computation  was  to  obtain  an 
approximate  bandwidth  of  the  fluctuations  of  the  radar  cross  sections.  The  cosine 
transform  technique  was  used  because  it  was  simple  to  adapt  to  machine  computation  and  it 
gives  the  closest  result  to  the  fast  fourier  transform  (FFT)  without  the  complexity  (12: 137- 
139).  Subroutine  COTRAN  performed  the  computations  to  produce  a  normalized 
magnitude  of  the  power  spectrum.  The  equations  used  in  the  cosine  transform  are 
presented  below. 


nr-i 

V  ...  . . 

spectr(K  At)  =  r  x(t)  cosvZJt  tk/nr)  K  =  u,l,...,nr/Z  (]/) 

i=l 

Af  =  l/(nr  *  T  )  (18) 

Next,  the  Rayleigh  and  Weibull  variables  were  generated.  The  derivations  of  these 
variables  were  discussed  in  chapter  two. 

At  this  point,  the  model  has  generated  a  set  of  data  and  periormed  its 
transformations  and  calculations.  The  last  major  section  of  the  program  dealt  with  testing 
the  variables  for  validation  purposes.  The  IMSL  subroutine  GITT  was  called  three  times  to 
test  the  goodness  of  fit  of  the  three  generated  vectors  y,  ray,  and  weib  to  their  theoretical 
gaussian,  Rayleigh,  and  Weibull  distributions.  A  graph  of  each  sample  vector  is  also 
plotted  vs.  the  theoretical  cumulative  distribution  function  with  95  percent  confidence  bands 
using  IMSL  subroutine  USPC.  The  IMSL  subroutine  VSRTA  merely  sorts  the  vector  into 
ascending  order  for  graphing  purposes. 
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The  program  also  converts  the  uncorrelated  vector  X  and  the  correlated  vector  X 
into  vectors  with  means  of  40.  This  was  done  to  demonstrate  a  chaff  cloud  simulating  a 
decoy  target  with  an  average  RCS  of  40  square  meters. 
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The  purpose  of  this  chapter  is  to  generate  several  sets  of  data  and  validate  the  model 
through  statistical  testing  and  by  making  observations  on  the  data. 

First,  the  program  generated  independent  gaussian  vectors  with  means  of  zero  and 
variances  of  one.  Three  tests  were  performed  for  each  value  of  nr  equal  to  20, 50,  and 
100.  Each  test  used  a  different  seeding  number  to  drive  the  IMSL  random  number 
generators.  Table  I  shows  the  results  of  these  tests. 


nr 

test 

20 

I 

999348 

.0000000 

2 

60259 

.0000000 

3 

565553 

.0000000 

50 

1 

39814055 

.0000000 

2 

24681055 

.0000000 

3 

872257 

.0000000 

100 

1 

7866543 

.0000000 

2 

629895 

.0000000 

3 

2833745 

.0000001 

1 

.961984 

.860671 

1.034006 


.967838 

.887949 

1.206430 


1.242569 

1.055229 

1.073428 


As  expected,  the  means  are  all  virtually  zero  due  to  the  fact  that  the  program 
subtracts  oui  the  mean  of  the  initially  generated  random  numbers.  However,  the  variances 
fluctuate  around  one  but  tend  to  approach  one  as  nr  increases.  This  foilows  from  basic 
probability  theory  that  as  the  number  of  samples  increases,  the  mean  tends  to  approach  a 
constant  number.  The  variance  also  tends  to  fluctuate  less  dramatically  as  the  number  of 
samples  increases.  So  the  model  behaves  well  in  producing  vectors  of  gaussian  random 
numbers  with  zero  mean  and  variances  of  one. 
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Next,  the  model  created  the  correlated  variables.  To  show  the  effect  of  the 
correlation,  two  graphs  were  generated  to  compare  the  independent  variables  with  the 


Time  (msec) 

Figuie?.  Correlated  vs.  Uncorrelated  Vectors  For  nr-20 
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Radar  Cross  Section  (m2) 


correlated  variables  for  nr  equal  to  20  and  50.  Figures  2  and  3  show  these  relationships.  A 
similar  effect  can  be  shown  for  nr  equal  to  100. 


The  effect  of  the  correlation  is  very  evident  and  thus  confirms  that  the  model  is 
indeed  smoothing  out  the  fluctuations  in  a  dramatic  way. 


0  10  20  30  40 


Time  (msec) 

Figure  3.  Correlated  vs.  Uncorrelated  Vectors  For  nr=50 
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The  next  section  of  the  program  computed  the  autocorrelation  of  the  correlated 
samples  to  provide  a  numerical  validation  of  the  data  graphically  represented  in  figures  2 
and  3.  Two  graphs  were  produced  for  nr  equal  to  20  and  100  in  figures  4  and  5. 


Figure  4.  Autocorrelation  vs.  Time  i^ag  for  nn=20 
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Figures.  Autocorrelation  vs.  Time  Lag  for  nr=  100 


K 


A  few  quick  calculations  can  show  that  using  equation  1  to  compute  the  variance  - 
covariance  matrix  will  generate  numbers  in  the  range  of  1  to  .90574  for  time  lags  of  0  to 
.099  seconds.  This  corresponds  to  99  time  intervals  :sing  a  1000  Hz  PRF  radar  as 
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specified  in  this  model.  Also,  the  numbers  generated  by  equation  2  will  decrease  from  1  to 
.999901.  Both  sets  of  numbers  decrease  in  a  nearly  linear  fashion.  This  may  seem  wrong 
at  first  since  both  functions  are  exponential.  But  it  is  easy  to  show  that  using  only  .099 
time  units  is  a  very  small  portion  of  the  exponential  curve  and  this  portion  near  zero  closely 
approximates  a  linear  function.  Thus,  the  autocorrelation  curves  in  figures  4  and  5 
accurately  portray  the  true  correlation  that  was  specified  in  the  beginning. 

The  graphs  also  show  that  the  vectors  exhibit  a  correlation  time,  the  time  required 
for  the  autocorrelation  to  decrease  to  zero,  that  is  related  to  the  number  of  samples  used. 
For  example,  19  msec  for  20  samples,  49  msec  for  50  samples,  etc.  This  is  due  to  the  fact 
that  a  limited  number  of  samples  was  "windowed"  out  of  an  infinitely  long  random 
process.  A  larger  number  of  samples  would  result  in  a  correlation  time  that  was 
independent  of  the  sample  length  as  expected.  These  values  are  still  reasonable  though 
taking  into  consideration  the  dynamics  of  a  chaff  cloud  and  the  random  nature  of  the  dipole 
orientation,  drift  velocity,  and  phase  variations.  Although  these  parameters  were  not 
directly  modeled  into  the  simulation,  the  results  show  a  realistic  behavior  in  the  relation  of 
one  cross  section  to  another. 

Furthermore,  a  comparison  was  made  between  the  autocorrelation  functions 
generated  by  equations  1  and  2  .  The  function  generated  with  equation  2  was  slightly 
greater  than  the  first,  but  the  difference  was  at  most  3-4%  and  represented  no  significant 
advantages  over  the  equation  1  function. 

The  power  spectral  density  (PSD)  was  computed  by  the  program  next.  Figures  6, 
7,  and  8  show  the  PSD  for  nr  equal  to  20,  50,  and  100  respectively. 


Again  it  can  be  seen  that  as  nr  increases,  the  graphs  tend  to  display  a  better 
representation.  In  this  model,  the  power  should  be  concentrated  around  zero  because  the 
PSD  represented  the  degree  of  fluctuations  that  the  mean  value  undergoes.  After 
correlation  of  the  random  variables,  the  fluctuations  should  be  very  small,  and  the  graphs 
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as  also  made  concerning  the  power  spectral  densities 
.lick  computation  can  show  that  the  values  generated 
than  the  equation  1  values.  This  is  not  surprising 


since  equation  2  is  just  the  square  of  equation  1.  On  a  normalized  magnitude  plot,  the 
values  are  less  than  or  equal  to  1  which  would  result  in  a  3  decibel  loss  when  the  value  is 
squared  instead  of  a  gain,  as  might  be  expected.  Again  no  significant  advantages  appear 
from  using  equation  2. 

The  last  section  of  the  algorithm  computed  the  goodness  of  fit  of  the  generated 
sample  gaussian,  Rayleigh,  and  Wcibull  random  variables.  Three  tests  were  again 
performed  with  different  seeding  numbers  for  each  nr  equal  to  20, 50,  and  100.  Table  H 
displays  these  results. 

TABLE  II 

Comparison  of  Goodness  of  Fit  Tests 


KTSTiutfHHi 

1 

| 

mmmim 

20 

56347856 

45678 

.6594 

.0786 

49 

65478 

.3618 

.4493 

.1116 

2938445 

864120 

.3080 

.8495 

.1991 

50 

705449 

8237126 

msm 

.0186 

.0719 

434819 

23569 

.0905 

.0992 

76239 

49876123 

.0160 

.0008 

100 

99267 

8839212 

.3286 

— 

.1847 

234598723 

102 

.1197 

■ 

.0215 

8092 

5454 

.5724 

■ 

.4481 

One  problem  arose  when  computing  these  values.  It  is  easily  seen  that  some  values 
are  quite  good  (approaching  0),  implying  high  accuracy  in  the  generated  variables.  But 
others  were  computed  around  0.5  and  even  higher.  It  appeared  that  the  chi-squared  statistic 
was  highly  dependent  on  the  seeding  value.  Also,  the  number  of  equiprobable  tally  cells 
influenced  the  chi-squared  statistic.  Cure  should  be  taken  in  specifying  this  parameter.  It 
was  observed  that  using  the  same  number  of  tally  cells  for  the  Weibull  and  Rayleigh  tests 
resulted  in  identical  chi-squared  statistics.  In  summary,  this  did  not  appear  to  be  a 
dependable  test  because  the  statistic  varied  so  much  with  the  input  parameters. 


An  alternative  test  performed  was  graphing  a  sample  vector  vs.  its  theoretical 
cumulative  distribution  function.  This  was  done  with  the  IMSL  subroutine  USPC.  This 
subroutine  plotted  the  sample  valuer  theoretical  values,  and  a  two-sided  95%  confidence 
band  around  the  sample  values.  Every  data  set  generated  produced  a  graph  that  was  well 
behaved  inside  the  confidence  bands.  These  graphs  are  presented  in  figures  9, 10,  and  1 1 
This  test  was  deemed  much  more  stable  and  suitable  for  confirming  the  distribution  of  the 
generated  variables. 

This  concludes  the  intended  validation  of  the  model  with  the  applicable  tests  and 
comparisons.  The  final  chapter  will  draw  some  conclusions  from  these  tables  and  figures 
and  make  some  suggestions  for  improving  the  model. 
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9.  Sample  Gaussian  Vector  vs.  Theoretical  Distribution  for  nr=100 


igure  U>  Sample  Wcibull  Vector  vs.  Theoretics’  Distribution  For  nr =100 


Figure  if.  Sample  Rayleigh  Vector  vs.  Theoretical  Distribution  For  nr=100 


V. 


Corclusions 

The  purpose  of  this  study  was  to  develop  a  computer  algorithm  capable  of 
generating  correlated  gaussian  random  sequences.  These  sequences  could  be  used  to 
represent  consecutive  chaff  cloud  radar  cross  sections  in  electronic  warfare  training 
scenarios.  Several  sets  of  data  were  generated  by  the  model  and  subjected  to  various  tests 
and  observations. 

The  model  proved  to  be  well  behaved  in  all  areas  tested  except  the  numerical 
goodness  of  fit  tests.  Several  inconsistencies  were  noted  involving  the  relationship 
between  the  input  parameters  to  the  IMSL  subroutine  GFTT  and  the  resulting  chi-squared 
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varied  so  dramatically  with  the  input  seeding  numbers  and  the  number  of  tally  cells 
specified.  This  problem  /as  offset  by  using  a  graphical  method  of  determining  the 
goodness  of  fit.  The  data  produced  by  the  model  was  declared  satisfactory  by  falling 
within  the  95%  confidence  bands  on  all  data  generated. 

The  model  realistically  generated  the  independent  gaussian  variables  required  and 
the  corresponding  Rayleigh  and  w'eibuii  variables  that  were  derived  from  the  gaussian 
numbers  all  proved  to  be  well  behaved.  The  large  fluctuations  present  in  the  independent 
variables  were  transformed  into  smoothly  varying  ones  just  as  the  model  was  designed  to 
do. 

Furthermore,  the  autocorrelation  and  power  spectral  density  properties  followed 
from  the  predicted  theories.  The  autocorrelation  should  have  been  linear,  and  the  PSD 
should  have  been  concentrated  around  zero.  It  was  noted  as  nr  increased,  these 
computations  became  steadier  and  more  well  behaved. 
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Ill  summary,  this  model  accomplished  all  of  its  initial  objectives  satisfactorily. 
Several  improvements  can  be  made  and  there  are  areas  which  could  be  studied  more 
thoroughly  as  discussed  in  the  next  section. 

RccomiDfiadations 

Further  study  should  be  done  into  the  problems  discussed  concerning  the  goodness 
of  fit  tests.  An  analysis  might  be  performed  to  determine  the  optimum  number  of  tally 
categories  that  would  result  in  the  most  consistent  chi-squared  statistic  results. 

A  study  could  be  done  to  see  if  an  upper  or  lower  limit  exists  on  the  PRF  that  may 
inject  some  inconsistencies  into  the  model. 

Determine  if  there  is  a  more  efficient  or  more  accurate  algorithm  available  for 
computing  the  power  spectral  density. 

Consider  what  additional  tests  the  model  could  be  exposed  to  that  might  simulate  a 
training  environment  and  determine  if  the  model  is  realistic  in  simulating  real  world 
dynamics. 

The  model  could  be  used  as  suggested  earlier  to  test  other  radar  systems  and  their 
effectiveness  against  chaff,  or  it  could  be  used  to  evaluate  a  model  of  a  radar  system. 

Introduce  the  blooming  effects  of  a  dispensed  chaff  cartridge  into  the  model. 
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Appendix  A  -  Probability  Theory 


In  this  appendix,  a  brief  review  of  probability  and  statistical  theory  will  be 
presented.  This  section  is  included  to  make  the  entire  paper  self  contained. 

A  random  process,  stochastic  process,  or  time  series  is  a  random  variable  that  is  a 
function  of  time  (14:3).  A  random  variable  X  and  its  value  at  time  tk  may  be  denoted  as 


Xk  =  X(tk) 


(A-l) 


Because  of  its  random  nature,  it  is  meaningless  to  talk  about  the  value  of  a  random 
variable  obtained  at  a  certain  time  or  values  observed  during  a  particular  time  interval. 
Instead,  the  idea  of  probability  provides  precise  definitions  of  certain  distributions  and 
averages  which  can  be  predicted  and  observed  with  some  level  of  confidence  (14:3).  The 
first  and  second  order  probability  distributions  can  be  defined  as 

x 

F(x,t)  =  Prob  [X(t)  <  x  ]  =  |  f(v)  dv 

-oa 

F(x1,x2:tl,t2  )  =  Prob  [X(tj)  <  xlf  X(t2)  £  x2  ] 

xi  x2 

f(u,z)  dz  du 

-OQ  -OQ 


(A-2) 


(A-3) 


The  probability  density  function  (pdf)  f(x)  and  the  cumulative  distribution  function 
(cdf)  F(x)  are  related  more  dearly  in  equations  A-4  and  A-5. 


a 

Fx(a)  =  Prob  (--»  <  x  <  a )  =  1  f(x)  dx 

-w 


(A-4) 
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f(x)  =  d/dx  [F(x)j 


( 

i 
j 

j 

(A  -5) 

Oftentimes,  the  pdf  of  a  random  variable  X  is  known  and  the  pdf  of  a  function  of  x, 

1 

for  example  y=f(x)  is  desired.  If  the  inverse  function  x=f'l(y)  is  single  valued,  then 
fy(y)  53  fx(x)/(dy/dx|  (A-6) 

The  formuia  holds  for  multiple  valued  functions  (13:95)  such  as 

-  f.-Ky),  x2  =  f^-’fy), ... 

fy(y)  -  fxiXil/ldy/dxJ  fx(x2)/|dy/dx2i  +  ...  (A-7) 

Random  variables  are  usually  classified  as  stationary  or  nonstationary  depending  on 
whether  the  statistics  of  the  variables  are  independent  of  or  dependent  on  the  independent 
variable  lime.  Since  this  paper  deals  strictiy  with  wide  sense  stationary  (w'SS)  variables, 
only  that  type  will  be  described  here.  A  variable  that  is  WSS  possesses  a  first  order  pdf 
that  >s  independent  of  ume  (a  constant),  and  its  second  order  pdf,  known  as  the 
autocorrelation  function,  is  dependent  only  on  the  time  difference  (14:4). 

The  expected  value  or  mean  of  a  random  variable  X  is  defined  as  follows. 

oo 

r 

E[  X  1  =  =  J  x  f(x)  dx  (A-8) 

-oo 

The  variance  and  the  standard  deviation  of  the  random  variable  X  ate  defined  as 

VAR(X)  =  s2  =  E(X2)  -  E-fX)  (A-9) 

SD(X)  =  IVAR(X)11/2  (A- 10) 

where  E[  X  3  is  commonly  referred  to  as  the  DC  component,  E(X2)  is  the  mean  square 
value,  and  SD(X)  is  the  AC  component  if  X  represents  a  voltage  ( 14:5-6). 
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One  of  the  most  important  second  order  quantities  of  interest  is  the  autocorrelation 
function  R(t),  which  is  defined  as 

R(t)  =  E(2Li  Ka)  =  E[x(t)  x(t+T)j  (A-ll) 

It  may  be  noted  that  for  t=0,  R(t)  is  equal  to  the  mean  square  value,  and  as  r 
approaches  infinity,  R(t)  approaches  the  DC  value  squared  (14:6). 

The  power  spectrum  or  power  spectral  density  (PSD)  is  another  second  order  pdf 
that  is  of  importance.  It  is  the  fourier  transform  of  the  autocorrelation  funcuon  and 
provides  information  about  where  the  energy  in  a  function  is  concentrated  (13:265).  It  is 
defined  as  follows. 


X(f)  =  j  R(r)  exp(-j2it  fir)  dr  (A-I2) 

-CO 

nr-1 

X(f)  =  t  ^  x(i)  exp(-j  2rcif r  )  -<»  <  f  <  »  (A-13) 

i-0 

It  can  easily  be  shown  that  since  R(-t)=R*(t),  the  power  spectrum  of  a  real  or 
complex  valued  process  is  a  real  function  of  f.  Other  transforms  are  used  to  provide 
information  about  the  PSD,  including  die  cosine  transform.  The  cosine  transform  uses  this 
real  value  property  to  simplify  the  representation.  It  deals  directly  with  the  discrete  sample 
values  of  the  random  variable  and  does  not  require  transformation  of  the  autocorrelation 
function  (12:139). 

The  moment  generating  function  (MGF),  as  presented  in  chapter  2,  is  useful  in 
transforming  a  random  variable  from  one  distribution  to  another.  For  a  variable  Z=X2,  the 
MGF  is  defined  as 
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Mz(t)  =  E[exp(tz)]  =  E[exp(tx2)] 


oo 

J  exp(tx2)  f(x)  dx 

-©o 


Mz(t)  =  (l/v'Sw  a)  J  exp(tx2)  exp[-(x-p)2  /2a2  ]  dx 

«oo 


(A- 14) 


(A- 15) 


for  a  gaussian  random  variable  X  with  mean  p  and  variance  a2. 

After  expanding  the  exponential,  completing  the  square  and  making  some 
substitutions. 


Mz(t)  .  (1  -  2a2  t)'1/2  -  [  (l/2a2)  /  (l/2a2  - 1)  ]1/2  (A-16) 


T  f  /t2_  1  Aak  V?.  Urwt  «  4  Jljtfa.'L..*'  _  ..  <.L  _ _ A  _  _ _ _  C  C _ _  J  _  _ 

v  WIW44  uaa  iaii-a4U4iiAA  uidiuyuuuu  wuu  vnc  ucgioc  Ui  liccuum. 

Meyer  proved  that  the  sum  of  n  chi-squared  variables  each  with  l  degree  of  freedom  is  a 
chi-squared  variable  with  n  degrees  of  freedom  (9:201).  It  can  also  be  shown  that  this 
variable  is  exponential  with  parameter  a=l/2.  Below  are  moment  generating  functions  of 
some  well  known  distributions  (9). 


Uniform  on  [a,'o] 

Mz(t)  =  [exp(bt)  -  exp(at)]  /  (b-a)t 

(A- 17) 

Binomial  (n,p) 

Mz(t)  =  (p  cxp(t)  +  (l-p)]n 

(A- 18) 

Exponential  (a) 

Mz(t)  =  a/(a-t)  t  <  a 

(A- 19) 

Gaussian  (p,  a2) 

Mztt)  =  exp[pt  +  a2  t2/2] 

(A- 20) 
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Appendix  B 


Program  Listing 
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SIMULATION  OF 


*«r*»n 


CHAFF  CLOUD  *• 

SIGNATURE  ** 

** 

Captain  Richard  P.  Fray,  USAF  *" 

** 

l»********W#WW*W*«*«***»*****K***W*7****W»** 

r****tt***«*****w«tt*«********««****«*«*%***i» 


rea  1  ac(  1  00) ,  c<  100 ,  1 0SS )  ,  cell  a  (  S  > ,  comp  (  5  >  ,ca,freq(51  )  ,  q 
*  ,  q  1  ,  q2  ,  ray  (  IJJjBT)  ,rhol  100)  .  a  (  100  >  ,  amean  .  spe(  5  i  )  ,  «pactr  <  5 
*1 ) .sum. sum  1 .svar .sl< 100 > . vc ( 100 , 100) ,w(4jMT> ,weib< 100) , 
+  x ( 100) , xmean , xvar , xl < 100) , y ( 100) , ymean , yvar , y 1 ( 100 ) , s , 
*z 1 , z2 . zmean 

1 ntegsr  Ic,  idf  , ler. 1p,k,nr,nlZ,n95 

double  precision  daeed 

external  norm . we  1 b I , r ay  1 

data  cell  3 , f r eq . rho . spectr , x/ 307 *0.0/ 

!c»I 

idf -0 

ip-l 

nl2-2 

n9S-95 

sum-0.0 

auml "0.0 

tau  is  the  time  interval 
consecutive  radar  retv 

tau» .001 

nr  <s  the  numh 
cross  sectlr 
4  7 

nr» 100 

zme' 

S'' 


10 

15 

1 
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yi<  lm)«y(  fm)izmean 
sum'sumtyl  Im) 

16  continue 

ymean-sum/nr 

*  vector  yl  now  contains  variables 

*  with  a  Mean>an*in 

*  generate  a  second  vector  of  Independent 

*  random  variables  for  use  In  generating 

*  the  weibull  and  raylelgh  variables 

dseed-283153 
wr ite(  6, IBldseed 
18  format!"  daeed  2  *“,fl2.1) 

call  ggnm I (dseed . nr . s 1 
sum-# .0 
suel '0.0 
do  20  I d'l. nr 
sue'sum+s ( I d  ) 

20  continue 

smean'sum/nr 
do  25  I f « l . nr 

suml  »suml  ♦  <  s<  If > -smaan ) **2 

25  continue 
svar 'suml /nr 
sum-0.0 

do  26  In'l.nr 
s(  lnl«s(  ln)-ssiean 
sl<  fn>“s<  ln)+ziaaan 
uum-sua+9 ( I n ) 

26  continue 
ssaan'ius/'  - 

*  vector  si  now  contains  variables 

*  with  a  mean-zmean 

*  generate  the  desired  var lance-covar I ance  matrix 

do  32  1  ■ 1 , nr 

do  30  J'l , I 

If  <  I z . eq . 2  )  then 

veil ,J)«axp(-(abs( 1-J) )'tau**Z> 
e  1  ss 

vc<  i,J)«exp<-<abs< l - J ) ) "tau  > 
end  I  f 

30  continue 

32  cant  1 nuc 

*  generate  transformation  matrix  C 

do  40  ik>l ,nr 
c (  I k , 1) »vc(  I  k  .  i  ) 

40  continue 

do  48  lm»2,nr 
do  46  jm»2 ,  1  m 
If  i  i m . eq . jm  )  then 
.sum»0.0 

do  42  km- 1 .  I m- 1 
3um-sum+(c( im.km) )**2 
42  cont  I  nut 

c  <  tm. iml»sqrt(vc(  im,  t  m ) -sum  J 
e  1  se 

sum  1-0.0 

do  44  kn- l , jm- l 

sgml=»suml+c(  Im , kn > "c< jm , kn ) 


40 


44  continue 

c  <  lm,jm)-(vc(  Im.Jml-sumll/cIJm.Jm) 
end  1? 

46  continue 

48  continue 

*  generate  cor.  slated  v»ctor  X 

do  52  11  *  1 , nr 
do  50  Jk»l  ,  f  1 

x< 11 )-x( 1 1 >  +  <c< 1 1 . Jk> )*y< Jk  > 

5j0  continue 

52  continue 

*  calculate  mean  of  correlated  vector  X 

sum-0 . B 
do  6#  1  a*  l  , nr 
sum-sum*x \  I  a  > 

60  continue 

xmean-sum/nr 

*  compute  variance  of  correlated  vector  X 

sum  1-0.0 
do  70  Id- 1  . nr 

suml-suml+lxl 1 d > -xmean  )  »*Z 
70  continue 

xvar-iuml /nr 

*  translate  X  into  a  random  variable 

*  V 1  th  «  Bijsejsisr; 

do  72  1 4- 1 . nr 
x 1 1! 4  I «x ( t 4 ) -xmean+zmean 
72  continue 

*  calculate  the  autocorrelation 

*  of  the  correlated  vector  X 

do  32  1-1  ,nr 
sum-0.0 
do  30  J- 1 , nr 
JJ- J-1-1 

sum-sum+K '  J  >  *x  <  J  j  > 

80  continue 

ac 1 1 ) *sum/nr 
82  continue 

wr 1 1#( b . 8  4  > xmean . xvar 

84  formstl "  xmean* " . f 10. 7 , “  xvar - " , f 10, 7 ) 

do  38  k l « 1 , nr 
r ho ( k 1  1  -ac ( k l ;/ac( i > 

*  print  out  correlated  vector  and 

*  autocor re  1  at  1  on  values 

wr lte<  6,86>kl .xlkl > .kl ,xl <kl > ,kl-I .rholkl > 

86  formate  x  (  "  .  1 3 .  "  >» '  , f  9 . 5  .  -  x  l  (  "  .  t  3  .  "  >  «  ’ ,  f  9 . 5  .  "  rho(  " 
*. f 3. ’ >»" .f=9.5> 

88  continue 

*  compute  the  power  spectral  density 

*  of  correlated  vector  X 

call  cotran<nr,tau,freq.spectr.x,spe> 
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*  print  out  autocorrelation  values 
do  92  1*1, nr/2+1 

wr 1 te( 6,90)1-1 .freql 11,1-1 ,s pel  1  ) 

90  formatl"  f  r«q  l  *  ,  1  3  .  *  )  -  ’ ,  f  9 . 3  .  *  H*  spe<  *  ,  1 3 ,  *  >-■ 

*,f IB. 7) 

92  continue 

If  (  t  z  .  eq .  2  )  then 
goto  1  7jBT 

*  recompute  parameters  with  different 

*  var lance-covar lance  matrix 

el  se 

*  generate  welbull  and  rayletgh  variables 

do  102  l *  1 , nr 

*■' y< 1 ) > "*2  +  <  all)) ""2 

ql*l .0/3,0 

Z 1 “Z"*q 1 

we  1 b (  I ) *x  1 

qZ-0.5 

z2*z**q2 

riy(  1  )  *>z2 

*  print  out  Independent  normal  vectors  y  and  s, 

*  print  out  their  cor respond  1 ng  mean*xmean  vectors, 

*  and  print  out  welbull  and  raylelgh  vectors 

wr lte(6. 100) 1 ,y< I > , I .yl< 11 , I ,s( I > -  1 .weibi 1 ) . 1 ,rayl 1 1 
100  format! -y< " , 13." >«",f9.5."  y  1  ( *  ,  1 3 , “ > - ■ . f 9 . 5 . •  si*. 
+13," )-*,f9.5,"  welbl",  I  3  ,  *  )«”  ,f  9  .  S,  •  r«y(M3.  "  )«■, 
-f9.5> 

102  continue 

wr 1  tot  6 , 104 ) ymean , yvar . smean ,  svar 
104  f ormat ( " ymean* ’ , f 10. 7 . “  yvar ■ " , f 10. 7 . "  smean*", 

♦f 10 . 7 , “  svar-" ,f 10. 7! 

*  perform  goodness  cf  fit  test  on  gauss  Ian  variables 

k«4 

call  g-ritln  or  m.k.y,  nr,  cells,  comp  ,cs,  Idf.q,  ler) 
wr  1  te<  6 .  1 14) 

114  format  1 " goodness  of  fit  test  for  gauasljn  variables") 
wr  I  to  <G.116;cells<l).cells(2>,  cel  1  S  *  3  1  ,  Cfi  1  1  S  (  4  , 

*ca 1 1 s \ 5  1 

1  16  format ( " cell s-1 " . f 5 . 1 . ” , “ . f 5 . I . " , " , f 5 . 1 . * , “ , 

+f  5  .  1  .  "  .  "  .  f  5 . 1  .  “  > "  > 
wr 1 te ( 6 , 1  i  8  )  cs 

118  format (" ch 1 -squared  statistic  for  gauss  Ian  variables 
♦-  “.f9.5> 
wr  I  te ( 6 . 1 20  > 1 df 

120  formatl "gauss  Ian  degrees  of  freedom  *"  ,12) 
wr 1  tel  6 , 1 22  ) q 

122  f ormat ( “ probab I  1  1 ty  of  chi-squared  statistic  encsedln 
*g  gauss  Ian  C3  if  null  hypothesis  Is  trus  *",fl0.7) 

*  graph  normal  vector  vs.  theoretical  cdf 
call  vjrtal y.nr I 

cal  1  uspc(norm.y,nr,nl2.n95, <p,  (c.w) 

*  perform  goodness  of  fit  test  on  welbull  variables 
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idf-e 

k-4 

do  124  (-1,6 
calls! 1 >-0.0 
124  conilnut 

call  gf 1 t!we  fbl ,k .walb . nr , cel  1 s , comp , cs . Idf  ,  q  ,  lar  ) 
wr 1ta<6, 130) cal  Is! 1 > .calls! 2) .cal ls<3> .cal ls<4> , 

+  ca 1 1 s ! 5  > 

130  format(*c«lls-(",f5.1,“,*,fS.l.*,".f5.1,“,", 

+f 5 . 1 , “ , " , f 5 . 1 ,  *  )  *  ) 
wr Its! 6 . 132)cs 

132  format! *ch ( -squarad  stat(st(c  for  waibull  varlablas  ~ 
♦  * , f9.S> 
wr(t#(6. 134) (df 

134  f ormat ! “wa Ibu 1 1  degraes  of  fraadom  »  ",12> 

wr 1 ta! 6 , 1 36  >q 

136  format! 'probata ! 1 ity  of  chi-squarad  statistic  axcaadln 
*g  waibull  cs  (f  null  hypothesis  Is  true  “  ",fl0.71 

*  graph  waibull  vector  vs.  theoretical  cdf 
call  vsr ta ! wa lb . nr ) 

call  uspc ( wa lb  1 .walb.nr.nl2.h96.lp. lc.w> 

-  perform  goodness  of  fit  on  raylalgh  variables 

1  df -0 
k-5 

do  138  1-1,6 
calls! 1 >-0.0 
136  continue 

call  gfit(rayl.k,ray,nr,cel Is, comp , cs . i df . q , lar  > 
wr (te<  6 . 1 44 )ce1 1 s! 1 > , ce 1 1 s! 2  > . cel  1 s ( 3  > , cel  1 s! 4  ) , 

+c«l 1 s<  5  > 

144  format! -cal  Is-! “ . f 5 . I , “ , " . f 5 . 1 . ‘ . “ . f 5 . 1 . * ,  “  . 

. 1 . "  .  "  .  f 5 . 1  .  "  )  *  ) 
wr 1 te<  6 . 1 46  /  cs 

146  format! “ch  1-squared  statistic  for  raylelgh  variables 
+-  " . f 9 . 5 ) 
wr 1 te<  6 , 1 48  >  1  df 

148  format! "raylslgh  degrees  of  freedom  -  "  ,  1  2  ) 
wr 1 te<  6 , 1 50 >q 

150  format!  "orobab  1 1  ity  of  cM-squared  statistic  exceedln 
+g  raylelgh  cs  If  null  hypothesis  Is  true  -  “,fl0.7> 

*  graon  raylelgh  vector  vs.  theoretical  cdf 
call  vsrtalray.nr  l 

cal  1  uspclrayl .ray. nr ,nl2.n95, Ip, 1q ,w> 
end  1  f 

170  continue 
er.d 


*  subroutine  cotran 

*  computation  of  cosine  transform  t PSD >  of  samples 

* 

A'**wwtt*w****w******«**r*ww**»*******w*rt*"***w*****»cr«* 


w 

<r 


subroutine  cotr  an ! nr , tan . f req , spectr ,  x  ,  spe  > 
reel  f req! 51 > ,p 1 , spe! 5 1 > . spectr ( 5 1 > , tau , temp , thet . 
+  x< 100> , yy , zz 
I ntager  1  , ka , nr 
p 1-3. 141592654 
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do  205  ksi-l. nr/2+1 
kb«ka- 1 
t emp«0. 9 
yy*kb/ ( nr *tau ) 
f  reql  ka)-yy«-  .0001 
do  2 09  1«1. nr 
thet»Z*p1«1*kb/nr 
ii*cos( that ) 
temp "temp +2 i**!  1  ) 

200  continue 

*  spactr  Is  tha*  raw  computed  power  spectrum 
speetr  <  ka )"temp*tau 

*  spa  Is  tha  normalized  magnitude  of  the  power  spectrum 

spe<ka)aabs(spectr<ka)/spectr{  1  )  ) 

205  continue 
return 
end 


subroutine  norm 

gaussian  cumulative  distribution  function 


subroutine  nornlx.pl > 
pl».5*erfc<-x*. 7071066 ) 
return 
and 


subroutine  welbl 

weibull  cumulative  distribution  function 

'Wwr»*******««r*«*w*«**««rsf******w*«r****«w***WW***'K 


* 


* 


subroutine  we  I b 1 1 x , y I  ) 

y!»l-exp(-(x**3>f2) 

return 

and 


subroutine  rayl 

raylelgh  cumulative  distribution  function 

*************************** ***W*Wtt**WwW*W 


subroutine  rayHx.y2> 

y2«l~expt-<x*"2>/2> 

return 

erid 
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